# Appendix B
# Santiago Lopez-Cariboni, 2022. "Political Regimes and Informal Social Insurance", Comparative Political Studies.

# ---- wd ----
wd <- 'replication_materials' 
setwd(wd)

# ---- packages ----
pkgs <- c("data.table", "foreign", 'panelView', 'fect', "lubridate", 
    "ggplot2", "dplyr", "mFilter", "texreg", "AER", "xts", "interflex", "gridExtra")
suppressPackageStartupMessages(sapply(pkgs,require,character.only=T))
# Plot the effects (use Johannes Karreth's function: https://github.com/jkarreth/JKmisc/blob/master/ggintfun.R)
# devtools::source_url("https://raw.githubusercontent.com/jkarreth/JKmisc/master/ggintfun.R")
source('code/functions/ggintfun.R')
# Function for time Lags
source('code/functions/lagpanel.R')


# ---- data ----
dt <- fread('data/dt_replication.csv')
dt <- dt[LDC==1, ]


# ---- Figure B1 ----
# Left panel: Electricity Access
out.ife.access <- fect(elec.access.total ~ democracy 
    + l.imports
    + l.exports
    + log(l.pop)
    + log(l.gdp.pcap)
    + l.elec.cons.pc
    + l.pop.density, data = dt, index = c("iso3c","year"), 
  force = "two-way", method= "ife", CV = TRUE, r = 4, se = TRUE, 
  nboots = 200, parallel = TRUE, na.rm=TRUE)
pdf("figures/OA/B1_left_access.pdf", , width = 8, height = 6)
plot(out.ife.access, 
  main = "Access to Electricity (% of population)", 
  ylab = "Estimated ATT", 
  xlab = "Time since Democratization", 
  cex.main = 0.8, cex.lab = 0.8, cex.axis = 0.8)
dev.off()

# Right panel: TD losses
out.ife.losses <- fect(tdlosses ~ democracy 
    + l.imports
    + l.exports
    + log(l.pop)
    + log(l.gdp.pcap)
    + l.elec.cons.pc
    + l.pop.density, data = dt, index = c("iso3c","year"), 
  force = "two-way", method= "ife", CV = TRUE, r = 4, se = TRUE, 
  nboots = 200, parallel = TRUE, na.rm=TRUE)
pdf("figures/OA/B1_right_losses.pdf", , width = 8, height = 6)
plot(out.ife.losses, 
  main = "T&D Losses (% power output)", 
  ylab = "Estimated ATT", 
  xlab = "Time since Democratization",
  ylim = c(-15,15),
  cex.main = 0.8, cex.lab = 0.8, cex.axis = 0.8)
dev.off()



# ---- Table 1 ----
dt$l.elec.access.total <- dt$l.elec.access.total/100
m <- lm(outgap.tdl ~ l.outgap.tdl 
    + l.outgap.gdp.hamilton * l.democracy
    + l.elec.access.total
    - l.elec.access.total
    + as.factor(iso3c)
    + as.factor(year)
     ,
    data=dt[,], na.action=na.omit)
summary(m,  diagnostics = TRUE)
nc4 <- length(unique(m$model[['as.factor(iso3c)']]))

m.access <- lm(outgap.tdl ~ l.outgap.tdl 
    + l.outgap.gdp.hamilton * l.democracy
    + l.elec.access.total
    + as.factor(iso3c)
    + as.factor(year)
     ,
    data=dt[,], na.action=na.omit)
summary(m.access,  diagnostics = TRUE)
nc4.access <- length(unique(m.access$model[['as.factor(iso3c)']]))


m.access.controls <- lm(outgap.tdl ~ l.outgap.tdl 
    + l.outgap.gdp.hamilton * l.democracy
    + l.elec.access.total
    + l.capacity
    + l.imports
    + l.exports
    + log(l.pop)
    + log(l.gdp.pcap)
    + l.elec.cons.pc
    + l.pop.density
    + as.factor(iso3c)
    + as.factor(year)
     ,
    data=dt, na.action=na.omit)
summary(m.access.controls,  diagnostics = TRUE)
nc4.access.controls <- length(unique(m.access.controls$model[['as.factor(iso3c)']]))

v.names <- c(
    "Losses Cycle$_{t-1}$",
    "GDP Output Gap$_{t-1}$",
    "Democracy$_{t-1}$",
    "Democracy$_{t-1}$ $\\times$ GDP Output Gap$_{t-1}$",
    "Electricity Access (\\% of population)$_{t-1}$",
    "State Capacity$_{t-1}$",
    "Imports$_{t-1}$",
    "Exports$_{t-1}$",
    "Population (log)$_{t-1}$",
    "Real GDP per capita (log)$_{t-1}$",
    "Electricity Consumption$_{t-1}$",
    "Population Density$_{t-1}$"
    )



table1 <- texreg(list(m, m.access, m.access.controls),
    file= "tables/OA/T1_access_control.tex",
    label="tab:access_control",
    custom.model.names = c('(1)', '(2)','(3)'),
    caption="Ciclicality of Electricity Losses in Autocracies and Democracies controlling for Electricity Access",
    dcolumn = TRUE,
    no.margin=FALSE,
    fontsize="scriptsize",
    single.row=FALSE,
    use.packages=FALSE,
    booktabs = TRUE,
    digits=2,
    float.pos="htbp",
    sideways=FALSE,
    omit.coef= "(year)|(iso3c)|(Intercept)",
    include.rsquared = FALSE,
    stars = c(0.01, 0.05, 0.1),
    custom.gof.rows = list("Year dummies" = c("Yes", "Yes", "Yes"), 
        "Country fixed-effects" = c("Yes", "Yes", "Yes"),
        "Number of countries"= c(nc4, nc4.access, nc4.access.controls)),
    # reorder.gof = c(1, 2, 5, 3, 4),
    custom.coef.names=v.names
    )
